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Abstract. Motivated by a recent study casting doubt on the correspondence 
between Regge calculus and general relativity in the continuum limit, we explore 
a mechanism by which the simplicial solutions can converge whilst the residual of 
the Regge equations evaluated on the continuum solutions does not. By directly 
constructing simplicial solutions for the Kasner cosmology we show that the 
oscillatory behaviour of the discrepancy between the Einstein and Regge solutions 
reconciles the apparent conflict between the results of Brewin and those of previous 
studies. We conclude that solutions of Regge calculus are, in general, expected 
to be second order accurate approximations to the corresponding continuum 
solutions. 



PACS numbers: 04.20.-q, 04.25.Dm 

1. Introduction 

Regge calculus Q is a discrete theory of gravity which replaces the smoothly curved 
spacetime of general relativity with a lattice. The curvature of the lattice spacctime 
is concentrated entirely on the two-dimensional hinges of the four-dimensional lattice 
cells. 

Regge calculus holds much promise for the numerical investigation of both 
classical and quantum gravity. Although the lattice approach appears well suited 
to numerical applications, progress in the field has been slow. Only recently have the 
first completely generic four-dimensional numerical simulations been performed j^, ||| . 

After the recent papers by Brewin Q and M. Miller ||, and despite the proven 
track record of the Regge calculus, a lively debate arose || as to whether or not 
solutions of the Regge equations would converge to solutions of the Einstein equations 
in some suitable limit. Neither Brewin or M. Miller directly computed solutions of 
the Regge equations. Instead, they took the somewhat easier approach of evaluating 
the Regge equations on an exact solution of Einstein's equations. They did this using 
sophisticated interpolation schemes, based on geodesies, to map a range of Einstein 
solutions onto simplicial lattices. The residual of the resultant Regge equations, 
calculated using these interpolated lattice edge lengths, was then examined in the limit 
of very fine lattice discretisations. Brewin observed that the residual scaled as 0(1) as 
the lattice was refined, and he inferred that solutions of the Regge equations, in generic 
spacetimes, would not converge to solutions of the Einstein equations in the limit of 
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fine discretisations. M. Miller chose a slightly different measure for the residual of the 
Regge equations, but when expressed in terms of Brewin's measure, M. Miller's results 
are consistent with those obtained by Brewin. Despite this M. Miller does draw sharply 
different conclusions, based upon further calculations associated with averages of the 
Regge equations, namely that solutions of the standard Regge equations will converge 
as the fourth order in the lattice spacing to solutions of the Einstein equations. Similar 
behaviour, in which the residual converges only over averages of the equations, has 
been observed by Sorkin |?j in the context of massless scalar fields on a 2-dimensional 
simplicial space. 

The situation became even more confused after the work of Gentle and W. Miller 
j| , who demonstrated explicit quadratic convergence of particular solutions of Regge 
calculus to a solution, the Kasner spacetime, of general relativity. This spacctimc 
was one of the test cases used by Brewin and M. Miller. This seems most odd - we 
appear to have a convergent set of solutions from an apparently non-convergent set of 
equations. How can this be? 

In this paper we explore a possible explanation for this behaviour, as proposed 
by Brewin Q|, which is consistent with all previous numerical studies. It is important 
to note that the direct solution of the Regge equations has always generated solutions 
which converge to the corresponding solutions of general relativity. Although the 
observed rate of convergence is dependent upon the particular lattice and symmetry 
restrictions imposed, when the lattice construction allows the full expression of the 
gravitational degrees of freedom, second order convergence of the lattice solutions to 
the continuum has always been observed ||] . 

We begin with a more complete description of the problem in section |^, and then 
discuss a possible explanation in section |^. Finally, in section we present numerical 
evidence to support our proposal. 



2. Investigating the convergence of Regge calculus 

There are various ways in which the convergence properties of a numerical technique 
can be explored. The most rigorous approach involves the direct comparison of 
the approximate equations and the full system, concentrating on the leading order 
discrepancies. The underlying assumption is that the solutions of the approximate 
equations will converge at the same rate as the approximate equations themselves 
converge to the original system. 

The Regge and Einstein equations are too complex for such a direct approach to 
be beneficial. Not only are the equations inherently complex non-linear systems, but 
before such an analysis can proceed one must choose a lattice upon which to express 
the Regge equations, together with a coordinate system in which to write out the 
Einstein equations. Once these decisions have been made, we are faced with a more 
problematic choice: how are the two sets of equations to be compared? In general 
there are more Regge equations (one per lattice edge) than Einstein equations (ten 
per spacetime event). 

There is no clear prescription for averaging the lattice equations to obtain the 
correct number of Einstein equations in the continuum limit. Brewin Q and M. Miller 
H have both developed their own schemes for averaging the Regge equations. A 
direct assault on the convergence properties of Regge calculus is therefore likely to 
be applicable only on a specific lattice, with a particular choice of averaging in the 
continuum limit. That is, the results would indicate that, in some averaged sense, the 
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Regge equations converge to the Einstein equations. However, the result would not 
constitute a proof valid for all possible lattice choices. 

Given these limitations, it is natural to ask if similar results can be obtained 
by directly solving the Regge equations and comparing the lattice solutions to 
corresponding Einstein spacetimes. Brewin and M. Miller have investigated the inverse 
problem: given a known (analytic) solution of general relativity, and interpolating that 
solution onto a lattice, what can be said about the convergence of the Regge equations? 

The approach taken by Brewin and M. Miller was to introduce some discretisation 
process in the smooth manifold on which the Einstein equations are defined. This 
provides a way to map the smooth metric solutions to a lattice, giving a set of 
lattice edge lengths Le derived from the continuum metric. Both authors used 
geodesic lengths to map continuum information onto the lattice spacetime. These 
new "continuum" lattice edge lengths do not in general satisfy the Regge equations, 
but we can evaluate the residual 



which is an indication of how well the interpolated Einstein solution satisfies the Regge 
equations. It is this residual, using an appropriately chosen norm, which both Brewin 
and M. Miller considered j|, [§]. 

Brewin observed that the residual of the simplicial equations remained roughly 
constant as the lattice was refined on a fixed region of spacetime. This result led 
Brewin to question the validity of Regge calculus, since any useful numerical scheme 
must converge to the underlying solution of the partial differential equations as the 
resolution is improved. M. Miller observed second order convergence of the residual, 
for any smooth metric, whether or not they were solutions of Einstein's equations. 
Had he used the same norm for the residual as used by Brewin we believe he too 
would have observed 0(1) convergence of the residuals. 

It is important to note that the experiments of Brewin and M. Miller do not 
directly evaluate the convergence of the numerical solutions. Rather, they investigate 
the convergence of the lattice equations to the Einstein equations. As the lattice is 
refined, it is reasonable to expect that the interpolated Einstein solution will satisfy the 
Regge equations increasingly accurately. This, however, is not what Brewin observed 
numerically. 

3. A possible explanation 

The observations of Brewin are particularly puzzling in light of the many previous 
applications of Regge calculus. Almost every numerical application of the method 
has displayed convergence towards the corresponding continuum solution, with 
most studies indicating that numerical Regge calculus is a second order accurate 
approximation to general relativity. 

The edge lengths measured in the Regge lattice (L) and the corresponding 
interpolated Einstein edges (L E , obtained by assigning geodesic lengths calculated 
in the continuum) are related as 



with previous numerical studies suggesting that p = 2 (see Gentle M for case studies 
and a general review). Throughout this paper we assume that S is a typical length 
scale in the lattice, and noting that the edges themselves are of this magnitude, we 



r = |R(L B ) 



(1) 



L = L E + 0(6 P+1 ) 



(2) 
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say that edges which satisfy equation (|^) are p-order accurate approximations to the 
continuum solution. 

Brewin [Q has proposed a mechanism whereby the solutions of the Regge 
equations converge to the corresponding Einstein solutions, while at the same time the 
residual of the Regge equations evaluated on the interpolated Einstein lattice edges 
does not converge. The key to Brewin's proposal is allowing the functional form of 
the error terms (the discrepancy between the Einstein and Regge solutions) to depend 
on the discretisation scale 8. This can be clearly seen from a toy model. 

Suppose £ is a second order differential operator, and y{x) is a solution of the 
equation 

A/ = o. (3) 

Brewin considers the function 

y s (x) = y(x) + 8 2 f(x/8) (4) 

for some arbitrary scalar 8 and an arbitrary, though bounded, function f(x). 
Considering the difference between the solutions we find that 

\y(x)-y s (x)\ =0{S 2 ), (5) 

indicating that the solutions differ only by "second order" terms. Noting that ys(x) 
is a solution of some other related equation C$, 

£sys = 0, (6) 

where Cs is also a second order differential operator (though different from C) , we find 
that the "residual" of ys{x) with respect to the original operator is 

\Cy s \ = 0{l). (7) 

This toy model embodies precisely the properties observed by Brewin - second order 
convergence of the solutions || , with no corresponding convergence observed in the 
Regge equations when they are evaluated on exact solutions interpolated from the 
continuum 0]. The discrepancy between the two solutions in this toy model is seen 
to be a wave-like disturbance with frequency proportional to 1/8. 

In the continuum limit it is reasonable to expect that the discrete Regge equations 
(or a weighted average over them) approach a system of differential equations. 
Moreover, consideration of the Einstein equations leads us to expect that the limiting 
form of the Regge equations is a set of second order non- linear equations. Furthermore, 
previous numerical experiments suggest that Regge calculus is a second order method; 
we expect p = 2. The toy model considered above leads us to expect that if the lattice 
solutions differ from the continuum solutions by terms with frequencies proportional 
to 1/8, it is not unreasonable that the residuals of the Regge equations remain roughly 
constant as the resolution is improved. 

This is precisely the behaviour observed by Brewin ; the solutions converge even 
though the residual of the equations do not. This explanation relies on the existence 
of high frequency, low amplitude waves in the simplicial solutions; a possibility not 
ruled out by any of the largely low-resolution applications of Regge calculus to date. 
In fact, one study hints at the existence of precisely this type of wave-like structure in 
the Regge solutions [||. In the next section we construct high resolution solutions of 
the Regge equations in order to gain insight into the fine-scale behaviour of the lattice 
solutions. 
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4. Numerical solution of the Regge equations 

In this section we solve the Regge equations for the vacuum Kasner cosmology using a 
(3 + l)-dimensional formulation of Regge calculus described elsewhere The initial 
value problem is solved and the lattice is evolved subject to simplicial lapse and shift 
conditions. 

The Kasner cosmology is a homogeneous anisotropic exact solution of the vacuum 
Einstein equations, and is the prototype for generic velocity dominated cosmological 
singularities. It also provides a model for epochs in the mixmaster cosmology between 
bounces. The metric may be written in the form 

ds 2 = -dt 2 + t pi dx 2 + t P2 dy 2 + t P3 dz 2 , (8) 

where the field equations reduce to a pair of algebraic constraints, 

Pi+P2+P3=pI+P 2 2 +pI = 1- (9) 
The Kasner exponents (pi,P2,P3) are a one-parameter set, and the solutions include 
the flat-spacetime case (0, 0, 1). 

The previous application of Regge calculus to the Kasner cosmology was viewed as 
a benchmark of the numerical technique, and as such the convergence of the solutions 
was an important issue. It was found that the lattice solution converged to the 
continuum as the second power of the typical lattice spacing This earlier study 
used a fixed lattice resolution, containing 128 vertices arranged in two offset cubic 
grids each consisting of 4 x 4 x 4 vertices. Reducing the overall scale of the lattice, 
combined with the T 3 topology of the Kasner solution, allowed the convergence rate 
to be estimated. 

In this paper we directly subdivide the lattice whilst fixing the size of the spatial 
region. This is a more general and demanding convergence test, as any long wavelength 
modes excluded in the previous study are now resolvable. Since the major issue 
is convergence we refine the lattice in only one dimension, but reduce the scale 
of the remaining spatial axes to avoid the creation of long skinny triangles. Such 
long and skinny triangular elements are known to cause instabilities in finite element 
calculations. 

This effectively "one-dimensional" model allows convergence analysis to be 
conducted whilst containing the computational scale of the problem. All calculations 
are performed on a lattice containing two offset cubic grids which are fully subdivided 
into triangles, tetrahedra and four-simpliccs. Each cubic grid contains 2™ x 4 x 4 
vertices (N = 2 n = 4, 8, . . . , 1024), with the x-axis having a fixed length of X = 10, 
and the remaining axes scaled to compensate; Y = Z = 40/2™. This ensures that the 
spatial edges surrounding a vertex are all of the same order of magnitude. The typical 
spatial resolution in each model is thus S rj 10/2". 

The full four-dimensional two-slice initial value problem is solved at t = 1, after 
which the lattice is evolved to t = 2 using the (3 + l)-dimensional Sorkin evolution 
scheme and "geodesic slicing" conditions; that is, with unit lapse and zero shift. The 
time step is chosen for each resolution to give a Courant factor of dt/l m 0.2; for 
N = 512 we chose dt — 3.9 x 10~ 3 . However, the main results of this paper were 
found to be largely insensitive to the choice of timestep, provided it satisfied the 
Courant condition. For full details of the evolution and initial value algorithms, see 
Gentle and W. Miller g] . The results of the Regge evolution are compared with the 
continuum solution by performing a least-squares fit of the function 

L(t) = Lot Pr (10) 
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to the time evolution of each class of the axis-aligned spatial edges, L x ,L y and L z . 
In agreement with previous studies, the lattice solution is found to remain spatially 
homogeneous throughout the evolution to within roundoff error. The typical standard 
deviation of the edges also remained below the limit of numerical accuracy. 

Figure [l] displays the fractional difference between the early evolution of the Regge 
cosmology and the exact solution for the axisymmetric case with p\ = p2 — 2/3, p% = 
— 1/3. Although the average magnitude of the discrepancy is small, there is a clear 
oscillation superimposed upon a general trend. The figure displays data for three 
lattices consisting of 128 x 4 x 4, 256 x 4 x 4 and 512 x 4 x 4 vertices. 

An estimate of the rate at which the Regge solution converges to the continuum is 
given in figure ||. After performing a least squares fit of equation ([!(]) to the simplicial 
solution in the region t = [1,2], we plot the fractional difference between the fitted 
power p r and the analytic value p c for each spatial axis. It is clear that this power 
law fit approaches the continuum solution as the second power of the lattice spacing. 
Figure || shows the scaled error in the least squares fit, giving an indication of how 
well the power law (|Io|) describes the Regge solution. Again, we find second order 
convergence. 

The data represented in figure |^ contains additional information. It is clear from 
figure [l] that there is a strong wave-like oscillation superimposed upon the general 
error curve. Fitting equation ( |To[ ) to such data will draw out the averaged solution, 
while the error in the fit will reflect the magnitude of any discrepancies from that 
average behaviour. Thus we see that both the oscillations and the underlying trend 
in figure [l] converge to zero as at least the second power of the grid spacing. 

In light of the discussion in the previous section, the frequency of the oscillations 
evident in figure [l] is also of interest. It is clear from the figure that the frequency 
increases with the lattice resolution. Figure [l] shows this in a more quantitative 
manner, displaying the Fourier transform of a small segment of the fractional error. It 
is apparent that the frequency of the waves varies linearly with the number of vertices, 
or alternatively, that the wavelength of the oscillations is proportional to 1/N. 

Together, these (3 + l)-dimensional convergence results show that the numerical 
Regge solution is a second order accurate approximation to the corresponding solution 
of the Einstein equations. Moreover, the leading order difference between the two 
solutions contains high frequency, bounded oscillations. We have shown that the 
frequency of these waves is proportional to N (or 1/5) , while their magnitude reduces 
as l/N 2 (S 2 ). This is precisely the variation postulated by Brewin Q, and strongly 
supports his explanation of the apparent non- convergence of the Regge equations. 

Finally, we note that the Kasner solution displays only temporal oscillations, 
with the spatial three-geometries remaining homogeneous to high accuracy throughout 
the evolution. This is likely a result of the high degree of symmetry of the Kasner 
spacetime, in which the constant time spatial hypersurfaces are flat. In a more general 
setting we would expect to see similar high frequency, bounded oscillations in both 
space and time. 

5. Conclusion 

We have carried out a rigorous convergence study of a particular solution to the Regge 
equations, and shown that the solutions do indeed converge to the corresponding 
Einstein solutions in the limit of very fine discretisations. 
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Following a suggestion of Brewin we examined the behaviour of the error terms, 
and found convincing evidence for the existence of high frequency, low amplitude 
oscillations superimposed upon the continuum solution. Brewin previously suggested 
that waves of precisely this form in the simplicial solution could explain the lack of 
convergence observed in the residual of the Regge equations when they are evaluated 
on the interpolated continuum solution. 

Together these results suggest that solutions of the Regge equations are generally 
second order accurate approximations to the corresponding Einstein spacetimes, with 
the discrepancy between the two solutions consisting of high frequency, low amplitude 
waves. Although these waves prevent the residual of the Regge equations converging 
to zero when evaluated on interpolated Einstein solutions g, they do not affect the 
overall second order accuracy of the simplicial solutions. 

We expect that all generic Regge simulations in vacuum will produce similar 
results. That is, all simulations will contain high frequency oscillations, linked to 
the inverse of the discretisation scale, which reduce in magnitude as the second 
power of that scale. These oscillations do not appear to induce instabilities in the 
numerical evolution of simplicial lattices. Indeed, previous studies have indicated that 
the oscillations gradually decay as the evolution proceeds Q. However, M. Miller |^| 
has reported instabilities arising in linearised Regge calculus on asymmetrical grids. 
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Figure 1. The fractional error SL/L of a typical edge in the Regge lattice as 
a function of time. Results are shown for spatial edges aligned with the x-axis, 

using lattice resolutions of 128 X 4 X 4 ( ), 256 X 4 X 4 ( ), and 512 X 4 X 4 

( — • — ) vertices. In all cases the total length of the lattice along the rr-axis is 
10.0 units. The magnitude of the observed oscillations are proportional to 1/-/V 2 , 
while the frequency increases with TV. 




Figure 2. A least squares fit of the function Lot Pr (equation |l0|) is performed on 
the edge lengths at different resolutions. We plot the fractional error in the fitted 
power p r as a function of the number of vertices. The fractional error in the fit 
along the x and y-axes (A-points; p rj 2/3) are indistinguishable; the z-axis fit 
(x -points; p fa —1/3) follow the power-law behaviour. The fractional error in the 
average Regge solution can thus be seen to vary as 1/N 2 . 
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Figure 3. The error in the least squares fit of equation ( |10[ ) to the mean edge 
lengths is shown as a function of the number of vertices. The results for the 
x and y-axes are again indistinguishable (A), with the z-axis showing the same 
trend (x). This is a measure of both the accuracy of the functional form given 
in equation and the amplitude of the waves shown in figure [lj We again 

find that the error measure converges to zero as 1/N 2 , indicating that the wave 
amplitude reduces as at least 1/N 2 . 
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Figure 4. The Fourier transform of the data in figure 1 is shown for lattices 

containing of 128 X 4 X 4 ( ), 256 X 4 X 4 ( ), and 512 X 4 x 4 (— • — ) 

vertices. The frequency of the wave-like oscillations apparent in figure hi are clearly 
proportional to N (or 1/(5). A quadratic drop-off in the power is also evident. 



